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ABSTRACT 


A MODEL TO INVESTIGATE THE INFLUENCE 
OF SUSPENDED SEDIMENT ON THE IIASS TRANSPORT 
OF A POLLUTANT IN OPEN CHANNEL FLOW 

Raymond Scott Chapman 
Old Dominion University, 1977 
Director: Dr. Chin Y. Kuo 

The environmental impact of the transport of pol- 
lutants in open channel flow has for many years been of 
interest due to the continuous introduction of heavy 
metals, pesticides, herbicides, and other foreign sub- 
stances into natural waterways. In order to fully 
understand chls transport process, it is necessary to 
examine the significance of its individual components. 

In the present study an explicit t^ro-dimensional finite 
difference model, designed to investigate the influence 
of suspended sediment on the pollutant transport process, 
IS presented. Specific attention is directed tov;ard 
examining the role of suspended sediment in 1) the tur- 
bulent vertical transport mechanism in a stratified 
flow, and 2) pollutant uptake due to sorption. Results 
presented indicate that suspended sediment plays a maior 
role in the pollutant transport process, and subsequently, 
any meaningful attempt to model the fate of a pollutant 
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in an alluvial channel must account for the presence 
of a suspended sediment concentration profile. Similarly, 
the vertical and longitudinal pollutant concentration 
distributions provided by the model may be utilized to 
improve upon the predictive capacities of existing water 
quality models. 
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CilAPTER I 


INTRODUCTIOiJ 

The environmental impact of the transport of pol- 
lutants in open channel flow has for many years been 
of interest due to the continuous introduction of heavy 
metals, pesticides, herbicides, and other foreign sub- 
stances into natural waterways. In order to fully 
understand this transport process it is necessary to 
examine the significance of its individual components. 
Since the initial investigations of Elder (1959) , much 
research has been directed toward modeling the influence 
of various factors such as bottom sediment, "dead zones," 
and bedload on the pollutant transport mechanism. 

Shih and Gloyna (1969) utilized a one-dimensional 
convective-dispersion model with a sorption function 
to examine the influence of bottom sediment on the mass 
transport of radionuclides in streams. Comparison of 
predicted results with observed data from flume experi- 
ments led to the conclusion that the effects of bottom 
sediment on the pollutant mass transport process could 
not be neglected in a model simulation. 

A similar model was adopted by Thackston and Schnelle 
(1970) to predict the effects of "dead zones" (pockets of 
little or no flow along stream banks) on the pollutant 
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transport process. To examxne the storage effects of 
"dead zones," methods for estimating values of volume 
and residence time parameters were investigated. Agree- 
ment between experimental and computed results suggest 
that the methods for predicting the effects of "dead 
zones" on longitudinal dispersion are valid. 

Using a one-dimensional stochastic model developed 
by Yang and Sayre (1971) , Shen and Cheong investigated 
the effects of bedload on the mass dispersion of a 
pollutant. Considering an instantaneous injection of 
contaminated sediment particles, they suggest: 

Contaminants such as herbicides, pesticides, 
radioisotopes can attach to sediment particles 
and move as bedioad in a stream, and the trans- 
port and dispersion of contaminants such as 
radioactive waste, can be affected by the dis- 
persion characteristics of the sediment particles. 

Similar ideas have been suggested with respect to 
the influence of suspended sediment on the pollutant 
mass transport process. In an assessment of the trans- 
port of radionuclides by streams (Sayre et al . , 1963) , 

It is concluded that: 

Available evidence indicates that sorption of 
waste materials from solution by stream sediments 
is the rule rather than the exception. Experi- 
ments have demonstrated that sorption of radio- 
active components from dilute solutions may 
exceed 90 percent, and that the concentration 
of radioactivity on the surface of sediment 
particles may become many thousands of times 
as great as the concentration in the surrounding 
liquid. 
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Reporting on the presence of DDT and Dreldrin xn 

rivers, Brexndenberg and Lxchtenberg (1963) concluded 

from results presented by Berck (1953) that: 

. . . chlorxnated hydrocarbons are adsorbed on 
the suspended solxds xn rxvers. Thus, the sxlts 
common to some rxvers may effectxvely remove 
pestxcides from water. 

In an experxmental xnvestxgatxon of sorptxon of 
pestxcides by clay minerals, Huang and Liao (1970) 
reported that DDT, Dxeldren and Heptachlor are readxly 
sorbed by clay mxnerals from aqueous solutions. As 
the suspended solxds settle out and accumulate as botrom 
sedxment, hxgh concentratxon of pollutant xn the sub- 
strate may result. Subsequently, Huang and Lxao suggest 
that: 


tinder certain condxtxons, part of the sorbed 
pestxcxdes can be desorbed and released xnto 
the water where rhey are maxntaxned by a 
dynamxc equxlibrxum system. Consequently, 
pestxcxdal desorptxon provxdes a contxnuous 
supply of toxxc materxal to water and creates 
many serxous water pollutxon problems. 

Rexmer and Krenkel (1974) xn an experxmental study 

of the uptake of Mercury by suspended sediment concluded 

that: 


Because xnorganxc mercury bxnds wxth sands, clays, 
and varxous organxcs, the contentxon that mercury 
pollutxon xs transported xn our v;aterways by 
sedxment xs supported. For example, the appear- 
ance of mercury contamxnatxon xn Kentucky Lake, 
whxch is over 100 mxles from mercury contamxnated 
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Pickwick Lake, the only known source of mercury 
in the system, can be explained on this basis. 

Thus, in an attempt to shed some light on this potential 
environmental hazard, an explicit two-dimensional finite 
difference model designed to investigate the influence 
of suspended sediment on the pollutant mass transport 
process is presented. Specific attention is directed 
toward examining the role of suspended sediment in 1) the 
turbulent vertical transport mechanism in a stratified 
flow, and 2) pollutant uptake due to sorption. The 
model developed in this study provides a quantitative 
description of vertical and longitudinal pollutant con- 
centration distributions which may be utilized to improve 
upon the predictive capabilities of existing water 
quality models. 
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CHAPTER II 
MODEL FORMULATION 


Governing Differential Equation 


Based on rhe principle of conservation of mass, the 


two-dimensional unsteady convective-dif fusion equation 
for a neutrally buoyant, conservative pollutant is of 


the form (Pritchard 1971, Sayre 1968, Harleman 1967): 


3C 

9 (Cu) 

3 (Cv) 

8 

' 8C 

— ^ 

+ 



e 

9t 

9x 

3y 

8x ' 

^ 8x 



( 2 . 1 ) 


where 

C *= C(x,y,t) = concentration of pollutant, 

u = u(y) = time averaged velocity in the 
X direction, 

V = V (y) = time averaged velocity in the 
y direction, 

~ ~ lc>ngitudinal diffusion coefficient, 

ey = Sy(y) = vertical diffusion coefficient, and 
r = sorption function. 

Assuming 

1) steady, uniform, incompressible flow, and 

2) v(y) - 0, 
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equation (2.1) becomes 


8C 

3G 

— + ^(y) — “ = e 

at 

3x 


f SC 

+ — 

1 Sy(y) 

3y 

i / 


(2.2) 


The coordinate system for equation (2.2) is oriented 
such that X defines the longitudinal direction along 
the centerline of the channel with positive x downstream. 
The positive y direction is defined vertically upward 
with y=0 corresponding to the channel bottom, and t 
representing time. Defining the following dimensionless 
parameters : 

X 

D 

_ Z. 

D 

u(y) 

U = 

U* 

6j,(y) 

Ex = (2,3) 

DU* 

ey(y) 

T? — 


D 
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where 

D = depth of flow, 

U* = shear velocxty = /t^Tp ~ /ghS'^ = / gDS 0 

g = acceleration of gravity, 2 

S_ = S„ = slope of the channel bottom — 

° ® ll.49R^’/ 

Tq = boundary shear stress, 

p = density of the fluid, 

R = hydraulic radius of the channel, 

n - Manning's roughness coefficient, and 

V = mean velocity over the depth of flow. 


equation ( 2 . 2 ) is written in the dimensionless form: 


3G 3C 3 dc\ 

t U — = + *— E„ — " + r* (2.4) 

3T 3X 8X^ 9Y \ / 


where r* is a dimensionless sorption function. 

Velocity Profile 

The classical representation of the velocity profile 
in. a turbulent shear flow is the von Karman-Prandtl 
logarithmic velocity distribution law (vanoni 1941) . 

When expressed in terms of the mean velocity over the 
depth of flow, the velocity law is written as 
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K = von Karman's constant. 

In nondimensional form, equatxon (2.5) becomes 

V 1 

U = — - + (In Y t 1) 

U* K (2.6) 

Realizing that equation (2.6) is undefined at Y=0, the 
velocity at the channel bottom is prescribed to be zero 
in accordance with a "no slip" condition. Subsequently, 
the velocity profile calculation starts at Y=.l. 

Graf (1971) suggests that in the presence of sus- 
pended sediment, von Karman’ s constant ic tends to decrease 
with increasing concentration, with a subsequent increase 
in the magnitude of the velocity profile. It is fully 
recognized that alternative velocity relationships have 

t 

been proposed to account for the presence of suspended 
sediment; however, the von Karman-Prandtl velocity law 
k=.4 is adopted due to its theoretical applicability in 
investigating the turbulent diffusion coefficients. 

Vertical Diffusion Coefficient 

Considering a neutrally buoyant pollutant, the 
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exchange coefficient for momentum is defined as 

T 


®m 


3u 

Pay 


(2,7) 


where 


T = shear stress. 

Assuming a linear variation in the shear stress with depth 
and noting that u|=Tq/p, equation (2.7) becomes 



ay 

Differentiating equation (2,5) with respect to y, 
vertical velocity gradient is expressed as 


( 2 . 8 ) 

the 


in - 

ay KY (2.9) 


Subsequently, substitution of equation (2.9) into equation 
(2.8) yields 

[ y\ 

©ni - tc U*y 1 

\ dJ (2.10) 

In nondimensional form, equation (2.10) is written as 


= icY(l - Y) 
m 


( 2 . 11 ) 



10 


where 





( 2 . 12 ) 


In applying equation (2.11) to represent vertical 

diffusion in the present model, the effects of density 

stratification due to the presence of a suspended sediment 

concentration profile must be taken into account. 

Leendertse (1975) suggests that; 

In fluids with vertically stable density gradients , 
each vertical motion of water mass has to work 
against buoyancy forces induced by the density 
gradient. If the available kinetic energy of 
turbulent motion is insufficient to overcome 
this stabilizing effect, turbulence is inhibited 
and suppressed. As a consequence, the process 
of vertical momentum exchange will be lower than 
under the neutral stability (vertically stable) 
condition . 

To model the relationship between the vertical exchange 
coefficient for momentum and the vertical diffusion 
coefficient, numerous empirical equations have been 
proposed of the form: 


ey = 


e^a + a 


Ri) 


-b 


(2.13) 


where by definition, the Richardson's Number 


Ri = 


~<3 


By 


!tl 


2 


(2.14) 
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and a and b are constants greater than zero« From 
equations (2.13) and (2.14), it is clear that for a 
highly stratified flow in which there is a large density 
gradient, the Richardson’s Number becomes large, and 
subsequently, the vertical diffusion process is suppressed. 
In the present study, the values of the constants a and 
b are chosen such that 


e 

-I = (1+3.33 Ri)"^/^ 
ra 

as reported by Munk and Anderson (19 48) . 
equations (2.11) and (2.15), the vertical 
coefficient is evaluated by 

icY(l - y) 

E = 

^ (1+3.33 Ri)^^/^ (2.16) 


(2.15) 

Combining 

diffusion 


Similar treatments of the vertical diffusion process in 
a stratified fluid are discussed by Harleman and Ippen 
(1967) and Jirka et al . (1975). Although equation (2.16) 
defines a vertical diffusion coefficient of zero at the 
free surface (1=1) , experimental evidence ( Jobson and 
Sayre 1970) suggests that in reality the value of the 
vertical diffusion coefficient at the free surface is 
greater than zero. Consequently, in the present study. 
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the value of the vertical diffusion coefficient at the 
free surface is prescribed to be one-half of the cal- 
culated vertical diffusion coefficient at Y=.9. 

Longitudinal Diffusion Coefficient 

Historically, investigations of longitudinal diffusion 
in open channel flow have been concerned with evaluating 
constant one-dimensional dispersion coefficients. The 
relationship 

E 

— ’ = constant, 

liU* (2.17) 

where 

E = longitudinal dispersion coefficient, 
is considered to be significant in relating the longitu- 
dinal dispersion coefficient to the flow condition 
(Sayre 1968) . However, Glover (1964) reports that 
experimental values of equation (2.17) vary from 5-24 
in laboratory flume experiments and from 15-800 in canals 
and natural streams. This discrepancy suggests that 
dispersion coefficients evaluated in laboratory studies 
may not be representative of the dispersion process in 
natural streams. Tennekes and Lumley (1972) state that 
"turbulence is not a feature of the fluid but of fluid 
flow," which implies that a laboratory experiment may not 
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accurately predict dispersion without accurately modeling 
the turbulent characteristics of the flow. 

In an attempt to formulate an applicable represen- 
tation for longitudinal diffusion as a function of depth, 
an alternative is to adopt what is called the "Four-Thirds 
Law" as discussed by Harleman (1966) , Sayre (1968) , and 
Blair (1976). In general, the "Four-Thirds Law" is 
written 

e^ = (2.18) 

where 

a = constant, 

G = energy dissipation rate per unit mass, and 
L = characteristic eddy length. 

Turner (1973) suggests that in a steady flow char- 
acterized by a logarithmic velocity profile, the energy 
dissiparion rate, G, can be evaluated as a function of 
the depth by 

9u 

G(y) = 

3y . (2.19) 

Substitution of equation (2.9) into equation (2.19) yields 

U*^ 

G(y) = 

KY 


( 2 . 20 ) 



14 


Presently, the constant a, and the characteristxc 
eddy length L must be determined experimentally for a 
given flow condition. If, however, a time averaged value 
of L is assumed constant, a new constant a', is defined 
such that 

a' = a (2.21) 


Substitution of equations (2.20) and (2.21) into 
equation (2.18) yields 




(icy) 


( 2 . 22 ) 


Thus, given that a prescribed value of occurs 

at a specified depth, the constant a' is uniquely deter- 
mined and a diffusion coefficient distribution is extra- 
polated over the depth of flow by means of equation (2.22) . 
In the present study, the prescribed value of e^^(y) is 
specified at middepth or, at'Y=.5. However, if experi- 
mental values of a and L are known, equation (2.22) is 
applied directly. In nondimensional form, equation (2.22) 
IS expressed as 


®x 


Cl 




(2.23) 
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Suspended Sediment. Profile 

The steady state one-dxmensional vertical turbulent 
diffusion equation for suspended sediment is of the form 
(Shen 1970) ; 


3Cs 



3 

3y 


3Cs \ 



(2.24) 


where 

- average particle fall velocity, 

Cs = suspended sediment concentration, and 

Dy = vertical turbulent diffusion coefficient 
for suspended sediment. 


Direct integration of equation (2.24) yields 

3Cs 


~V„Cs = D 
T y 


ay 


(2.25) 


which states that for an equilibrium sediment concentration 
profile, the quantity of suspended sediment settling must 
be balanced by the turbulent diffusion upward in the 
direction of decreasing concentration. Expressing the 
vertical turbulent diffusion coefficient as a function 
of the momentum transfer coefficient (Jobson and Sayre 


1970) , 


Dy = ^ 


(2.26) 
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where the constant of proportionality is generally 
considered to be less than or equal to one. Graf (1971) 
concludes that for a fine particle size, 0 is equal to 
one, and for a coarse particle size $ is less than one. 
However in the present model, 3 is assumed equal to one. 
Substitution of equation (2.10) into equation (2.25), 
with the subsequent integration yields , 


Cs = 



(2.27) 


where Co is a constant of integration. Defining Cs=Ca 
at a reference depth y=a where, in the present model, 
a/D=.l equation (2.27) becomes 


Cs 

Ca 


a 

’i - y ) 


DY 

[1 - ~l 


am 

V 



(2.28) 


where 

Ca = reference concentration, and 



tcU* 

} 


Noting that the suspended sediment concentration given by 
equation (2.28) is undefined at the channel bottom or 
Y=0, and zero at the free surface or Y=l, the following 
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boundary conditions are imposed (Kuo, 1976) : 

1) the suspended sediment concentration at the 
channel bottom (Y=0) is equal to the reference 
concentration Ca, 

2) the suspended sediment concentration at the 
surface (Y=l) is linearly extrapolated such 
that for a mean particle size <.0001 cm the 
concentration is prescribed to be that of the 
adjacent lower layer, and for a mean particle 
size >.G3 cm the concentration is set equal 
to zero . 

Consequently, calculation of the suspended sediment 
concentration is carried out between Y=.l, and Y=.9. 


Sorption Function 

In studies of equilibrium sorption (Poinke and Chester 
1972, Huang and Liao 1970, Reimer and Krenkel 1974, Boucher 

t 

and Lee 1972) , the empirically derived Fruendlich equation 
IS found to be representative of experimental results. 

As applied to the sorption of pollutants onto suspended 
sediment, the Fruendlich equation is written 


W 


M 




(2.29) 


where 

W = weight of sorbed pollutant by M grams 
of suspended sediment, 

C = concentration of pollutant in equilibrium 
with the suspended sediment, and 


k,in = experimentally determined constants. 
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Solving equation (2.29) for the weight of sorbed pollutant 
and multiplying both sides by the suspended sediment con- 
centration , yields 

Ce = k Cs (2,30) 


where 

Ce - equilibrium concentration of sorbed pollutant. 
Denoting r to be the time rate of change of sorbed pollutant, 

8Cas 

-r = = -k (Cas - Ce) 

9t ® (2.31) 


where 

Cas = an average concentration of sorbed pollutant, 
and 

kg = mass transfer coefficient. 

Therefore, substitution of equation (2.30) into equation 
(2.31) yields 


3Cas 

-r = = k_(kCs - cas) 

3t 


which when nondimensionalized becomes 


~r* 


3Cas 

= K (kCs C^/^ - Cas) 
s 


(2.32) 


3T 


(2.33) 
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where 

K 

® Hi (2.34) 

Model Suininary 

In summary, the equations which make up the present 
model are as follows; 

3c ac a^c a / ac ' 

— t u~“ = + ““ Ey — + r* 

3T ax ax^ 3Y \ 


where 

aCas 

r* - - = _K (kCsCvr^ - Cas) 

3T ^ 


V 1 

U = — + — (In Y + x) 

U* (c 

a* 

E = ■ ' 

^ (kY)1/^ (2.35) 

s 

kY(1 - Y) 

g « ■ ■ 

^ (1+3.33 Ri) 
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Implicit solution of the two unknowns C, the ambient 
concentration of pollutant, and Cas, the sorbed concen- 
tration of pollutant requires the reduction of large 
systems of nonlinear equations; subsequently, in the 
present model equations (2.4) and (2.33) are uncoupled 

I 

by solving for Cas with C evaluated at the previous 
time step. Thus, the pollutant concentration field is 
explicitly computed subject to the following boundary 
and initial conditions : 


1) The ambient and sorbed pollutant concentration 
fields are initially set equal to zero, or 


C(X,Y,0) 0, and 

Cas(X,Y,0) = 0. 


X > 0 


2) A line source of strength is prescribed over 
the entire depth of flow at X=0, or 

C(0,Y,T) = C^. T > 0 

3) Mass transfer across the free surface is specified 
to be zero, or 



0 . 


Y = 1 


4) To simplify treatment of the effects of bottom 
sediment and bedload, a flux is prescribed to 
the channel bottom. Specifically, the pollutant 
concentration at the channel bottom is calculated 
as a linear function of the concentration gradient 
between Y=0 and Y=.l. 
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CliAPTER III 
COMPUTATIONAL METHODS 

Mass Transport Equation 

To simplify programming procedures and minimize 
computation time, an explicit finite difference scheme 
is employed. The difference approximations are a modified 
upwind differencing method, where the derivatives are 
approximated forward in time, backward in the first 
derivative in space, and centered in the second deriva- 
tives in space. Therefore, the truncation error of 
this method is of order (Ax), (At), (Ay)^. Overall the 
scheme is first order accurate, however, it possesses 
both the conservative and transportive properties. 

Roache (1972) in a discussion on the use of a first order 
accurate upwind method suggests that; 

. . . it is also possible to more accurately 
represent -a derivative by using a nonconservative 
method, but the whole system is not more accurate 
if one's criteria for accuracy include the con- 
servative property. 

Realizing that conservation of mass must be maintained 
in a pollutant transport model, the conservative property 
IS indeed an important criteria when assessing the accuracy 
of a numerical scheme. Defining that a finite difference 
method is transportive "if the effect of a pertubation 
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IS advecte<3 in the direction of the velocity," Roache 
goes on to say. 

Innocuous and obvious as this definition may 
read, the fact is that the most frequently 
used methods do not possess this property. 

All methods which use center space derivatives 
for the advection term do not possess this 
property. 

Furthermore Roache suggests that: 

The transportive property appears to be as 
fundamentally important, as physically signi- 
ficant as the conservative property. At least 
in the sense, upwind differencing schemes which 
possess the transportive property are more 
accurate than schemes with space-centered 
first derivatives . 

Thus, although mathematically first order accurate, upwind 
differencing schemes may be preferable in pollutant 
transport modeling. However, Roache also points out 
that by adopting a non-centered difference operator 
to approximate the advection term, an artificial or 
numerical diffusion coefficient is introduced of the 
form: 


ae(y) = ^UAX(1 - c) ' 


( 3 . 1 ) 


where 


ae = numerical diffusion coefficient, 

c = courant number = ^ 

AX 


AX = longitudinal space increment, and 


AT 


time increment. 
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Therefore, to account for the effects of the numerical 
diffusion coefficient, the value of ae is computed and 
subsequently subtracted from the longitudinal diffusion 
coefficient at every iteration. To verify this correction 
method, the linear one-dimensional convective-diffusion 
equation, 


ac 

ac 

a^c 


— ^ 

V — = 

E — 


at 

3x 

ax^ 

(3.2) 


was solved numerically and compared to a known analytic 
solution subject to the following boundary and initial 
conditions: 

C(0,t) = 1 0 

C(x,0) =0 X > 0 

C{«,t) =0 t > 0 

Dailey and Harleman (1966) present an analytic solution 
to equation (3.2) which is of the form: 

C [vx\ /x + Vt\ lx - Vt^ 

~ s: Jjexp — erfc + h erfc 

Ci \E / 2>/Et \2/Et I (3.4) 

where the complimentary error function, erfc (s) =l-erf (s) 
and the error function is defined as 

2 

erf(s) = — 1 _-2 

/n / e ^ dC (3.5) 



In evaluating the analytic solution presented in equation 
(3.4 ) , the error function is approximated by means of 
the following series (Abramowitz and Stegun 1964) : 
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erf(s) = 1 - (aiq + 

j + a^jg"* + asq^)exp(-s^) (3.6) 

where 

ai == .254829592 

= -.284496736 
= 1.421413741 
a, =-1.453152027 
as = 1,061405429 
p = .3275911, and 

q = 1/(1 + ps) 

Figure (3.1) illustrates the close agreement between the 
numerical and analytic solution of equation (3.2). 

Denoting the subscripts to represent the 

position of X, Y, and T, respectively, the finite dif- 
ference approximations used in the present model are as 
follows: 

3C 1 


8T 


AT 


(3.7) 
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U- 


3C Uj 


8X 


(C - C ) 

I/J,N 


(3.8) 


(E„ - ae)j 
(C 


X. 


ax‘ 


(AX) 




+ _ ) (3.9) 

I-l / J,N 


a ac 

- Ey— 

3Y \ 3Y 


2 (AY) 




%J-i) 


C ^ ) 

I / J* X /N 


(3.10) 


acas 

3T 


— (Cas -' Cas ) 

I,J,N+1 I,J,n' 


(3.11) 


In general, numerical stability is maintained provided 
that the time increment 


AT < 


2 (E^-ae) 
(AX) 


2E, 


(AY) 


+ 


u 


AX 


(3.12) 


where maximum values of the turbulent diffusion coefficients, 
artificial diffusion coefficient and longitudinal velocity 


are used 
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Rxchardson* s Number 

Investigation of the effects of a suspended sediment 
concentration profile on the vertical diffusion mechanism 
requires that the Richardson's Number be evaluated at 
each depth considered. Noting that the specific weight 
Y=pg/ equation (2.14) is written, 

9Y 

Ri = 2 

f 3u \ 

(3.13) 


Substitution of equation 


(2.9) into equation (3.13) yields 


9Y 

9y 



(3.14) 


Assuming that the volume of suspended sediment is negligible 
with respect to the volunxe water, the specific weight y 
is approximated as a function of depth by 


Y(y) 


Tw 

1 - Cs (y) X 10“® 


(3.15) 


where 

Yw " specific weight of water. 

Therefore, by adopting a second-order accurate difference 
approximation, the vertical density gradient is evaluated 
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9Y 

3y 


2Ay 


(3.16) 


Fxgure (3.2) is a plot of Richardson's Number versus 
depth illustrating the sensitivity of the computed 
density gradient to the differencing increment in the 
Y direction. The change in the Richardson's Number with 
mesh refinement is primarily due to the specification 
of boundary conditions at the channel bottom and free 
surface. Nonetheless, close agreement is maintained 
between Y=.2, and Y=.8 which is essentially the region 
of interest in the present study. 
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CHAPTER IV 
MODEL SIMULATION 

Numerical simulation is approached in three phases, 
namely 

1) Dispersion experiments 

2) Stratification experiments, and 

3) Sorption experiments. 

In all the above simulations, the following model para- 
meters are held constant: 

D = 20 ft 
V = . 5 ft/ sec 
At = 5 sec 
Ax = 50 ft 
Ay = 2 ft 

Figure (4.1) illustrates the extent of the model simulation 

where , 

n = Manning's roughness coefficient, 

Emag = value of the longitudinal diffusion 

coefficient prescribed a.t middepth, and 

ds - mean particle size. 

The first phase of simulation examines the response of 
the pollutant transport process to variations in flow 
condition without imposing a suspended sediment concen- 
tration profile. The flow condition is prescribed by 



29 


choice of a channel roughness as defined by Manning’s 
roughness coefficient, and the magnitude of the longi- 
tudinal diffusion coefficient specified at middepth. 
Noting figure (4.1), Manning's roughness coefficients 
of .01, .03, and .06 are adopted as typical values for 
natural streams (Chow 1959), with a range for the longi- 
tudinal diffusion coefficient of 10, 50, and 100 ft^/sec. 
In the present model, with the depth and mean velocity 
held constant, increases in the channel roughness results 
in increasing channel slopes, and subsequently greater 
shear velocities. Figure (4.2) is a nondimensional 
plot of velocity versus depth illustrating the increase 
in shear velocity with an increase in channel roughness. 
In addition, from equation (2.9) it is readily seen 

that with an increase in shear velocity, a subsequent 

\ 

increase in the velocity gradient occurs. A similar 
result IS also indicated when plotting nondimensional 
longitudinal diffusxon versus depth (Figure 4.3). 

The second phase of experimentation examines the 
effect of stratification, due to the presence of a 
suspended sediment concentration profile, on the tur- 
bulent vertical mass transfer mechanism. Noting figure 
(4.1) the flow condition for this particular investi- 
gation IS prescribed by choosing a typical Manning's 
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roughness coefficient of .03, and a magnitude of 50 ft^sec 
for the longitudinal diffusion coefficient. In order to 
vary the shape and magnitude of the suspended sediment 
concentration profile, reference suspended sediment 
concentrations of 100,, 500, and 1000 ppm and mean particle 
sizes of .001, .005, and .01 mm are adopted. Figures 

(4.4) and (4.5) illustrate the change in the suspended 
sediment concentration profiles with varying reference 
suspended sediment concentrations and mean particle sizes. 
Thus, the variation in the suspended sediment concen- 
tration gradient, results in notable changes in the 
Richardson’ s Number distribution with depth (Figures 
(4.6) and (4.7)). Comparison of figures (4.4) and (4.6) 
reveals that the magnitude of the Richardson’s Number 

IS highly dependent on the suspended sediment concen- 
tration gradient? however, the distribution of the 
Richardson’s Number over the depth of flow appears to be 
more sensitive to the magnitude of the velocity gradient. 
Similar results are indicated when comparing figures 

(4.5) , and (4.7)? however, for ds=.01 mm the decrease 
in suspended sediment concentration is so rapid that 
as the velocity gradient approaches zero no apparent 
increase in the Richardson’s Number occurs. Recalling 
equation (2.16), the greater the Richardson's Number 
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the more significant the reduction in the turbulent 
vertical diffusion coefficient. This result becomes 
obvious when figures (4.6) and (4.7) are compared with 
figures (4.8) and (4.9), respectively. 

The third phase of experimentation investigates 
the uptake of pollutant by suspended sediment due to 
sorption. The flow condition is the same as in phase 2, 
and in addition, the analysis is further restricted by 
selecting a typical reference suspended sediment con- 
centration of 500 ppm, and a mean particle size of 
.005 mm (Figure 4.1). The sorption function presented 
in Chapter 11 is expressed as a first order reaction 
in terms of an equilibriiim uptake of pollutant for a 
given suspended sediment concentration. Experimental 
results relating the equilibrium uptake of pollutant 
to suspended sediment concentration are well documented 
in the literature from which the Freundlich constants 
of the pollutants considered in the present model are 
obtained (Table 4.1). Sorption rates are expressed as 
rate constants derived by simply taking the reciprocal 
of an average time in seconds for equilibrium uptake 
of a given pollutant by clays, and or sand. To account 
for the heterogeneous composition of suspended sediment, 
the rate constants are treated parametrically and varied 
over three orders of magnitude. 



32 


TABLE 4.1 


Approximate Freundlich 

Constants 

Pollutant 

k 

1/iti 

Reference 

Mercury 

5 X lO”* 

1 

Rieraer & Krenkel 
(1974) 

DDT 

3 X lO"* 

3 

Huang & Liao 
(1970) 

Heptachlor 

5 X 10"® 

5 

Huang & Liao 
(1970) 
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CHAPTER V 

RESULTS AND DISCUSSION 

Results of the first phase of experimentation indicate 
that the model simulation is sensitive to the choice 
of both the channel roughness and the magnitude of the 
longitudinal diffusion coefficient. Figure (5.1) is a 
normalized plot of the pollutant concentration versus 
depth illustrating the influence of varying channel 
roughness on the vertical concentration profile. The 
increase in magnitude of the vertical concentration 
gradient with channel roughness is a direct result of 
the increase in the velocity gradient due to greater 
shear velocities. Normalized plots of the pollutant 
concentration versus longitudinal distance, figures 
(5/2a,b,c), are presented to illustrate the effect of 
an increase in magnitude of the longitudinal diffusion 
coefficient on the longitudinal distribution of pollutant. 
Two important features of the dispersion simulation 
are displayed in these figures. The first is that 
early in the simulation, or at t=200 sec, longitudinal 
diffusion initially dominates the dispersion process as 
indicated by the increase in magnitude of the concen- 
tration wave front with greater values of Emag. However, 
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as the simulation proceeds in time, the longitudinal 
concentration gradient decreases and advection dominates 
longitudinal diffusion as the system approaches steady 
state (Figures (5.2b,c)). 

Results of the second phase of investigation reveal 
that damping of turbulent vertical mixing due to strati- 
fication caused by presence of a suspended sediment con- 
centration profile, IS a prominent feature of the pol- 
lutant transport process. Figure (5.3) , a normalized 
plot of pollutant concentration versus depth, exhibits 
a pronounced vertical concentration gradient due to 
the lack of vertical mixing. Although the influence 
of a suspended sediment concentration profile is notable, 
figures (5.4a,b) illustrate the lack of significance 
of the variations in magnitude of the suspended sediment 
concentration and mean particle size on the degree of 
reduction of turbulent vertical diffusion coefficient 
corrected for the influence of stratification is less 
than 30% of the maximum momentum exchange coefficient, 
figures (4,8) and (4.9), it is apparent that irrespective 
of the magnitude of the suspended sediment concentration 
or mean particle sizes considered the damping of the 
vertical turbulent mixing mechanism is significant. 

Analysis of the sorption experiments indicate that 
the uptake of Mercury, DDT, and Heptachlor by suspended 



35 


sediment must be considered a significant factor in the 
pollutant mass transport process.. Mercury, the only 
inorganic compound investigated is affected least by 
the sorption mechanism; however, the uptake of mercury 
by suspended sediment 'is quite notable. Figure (5.5) 
illustrates the decrease in the magnitude of the vertical 
concentration profile with increases in the order of 
magnitude of the sorption rate parameter kg. Further- 
more, an important trend is seen in the decrease of 
pollutant concentration with depth due to increasing 
suspended sediment concentration. This feature is 
further illustrated in figures (5.6a,b), where normalized 
plots of concentration versus longitudinal distance are 
shown at Y=.2 and .5, respectively. These results sug- 
gest that the distribution of the suspended sediment 
concentration profile must be accurately represented 
when modeling the uptake of a pollutant by suspended 
sediment in natural streams. DDT and Heptachlor are 
organic pesticides which, by their chemical nature, 
are highly receptive to uptake due to sorption by sus- 
pended sediment. Figures (5.7) and (5.8) are normalized 
plots of concentration versus depth and longitudinal 
distance, respectively, illustrating the increase in 
uptake of DDT with increases in the order of magnitude 
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of the sorption parameter kg. Similar results exhibited 
in figures (5.9) and {5.10} for Heptachlor suggest that 
the uptake of organic pesticides by suspended sediment 
is indeed significant. 
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CHAPTER VI 

SUMMARY AND CONCLUSIONS 

In summary, a two-dimensional convective-dif fusion 
model has been applied successfully to investigate the 
effects of suspended sediment on the mass transport of 
pollutants in open channel flow. The analysis is 
restricted to examining the transport of a conservative, 
neutrally buoyant pollutant in a uniform turbulent 
shear flow. The numerical simulations performed lead 
to the following conclusions: 

A. Application of mathematical models to investigate 
the pollutant transport process in natural streams 
requires the judicious choice of channel roughness, 
and a reasonable representation of longitudinal dif- 
fusion. 

B. The suppression of the turbulent vertical mixing 
process, as a result of vertical density stratification 
due to the presence of a suspended sediment profile is 
shown to be significant. 

C. The reduction of pollutant concentration due to 
uptake by suspended sediment suggests that the sorption 
mechanism is an important component of the mass trans- 
port process. 
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The present study has clearly illustrated the signifxcance 
of a suspended sediment concentration profile on the 
pollutant mass transport process. In order to refine 
the results presented herein, further research is needed 
to improve upon the following: 

A. As discussed in Chapter II, a von Karman's constant 
ic==.4 IS adopted without correction for the presence of 
suspended sediment. However, if the suspended sediment 
size distribution is known, correction can be made by 
adjusting von Karman's constant according to the 
empirical relation presented by Einstein and Adbel-Aal 
( 1972 ) . 

B. A theoretical approach to extrapolate a longitudinal 
diffusion coefficient profile has been proposed; however, 
experiment evidence is needed to test the applicability 

of the "Four-Thirds Law" to turbulent diffusion in natural 
streams . 

C. Constant sorption rates are adopted in the present 
model, where in reality sorption rates are decreasing 
nonlinear functions. A possible means to more accurately 
represent sorption rates would be to develop an empirical 
rate function from experimental data. 

D. In the present study, the analysis has been restricted 
to examining the fate of conservative pollutants. To 
account for the decay of a nonconservative pollutant 
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such as radionuclides, it is only necessary to subtract 
the concentration of pollutant computed by an experimental 
decay rate function from the values of C and Cas at every 
time Iteration. 

E. A constant longitudinal suspended sediment concen- 
tration at each depth is assumed throughout the present 
model. However, in applying the model to a physical 
situation, remote sensing techniques may be employed 

to provide a longitudinal distribution of surface sus- 
pended sediment concentrations with a minimal amount of 
field data for ground truth. Thus, by redefining the 
suspended sediment concentration reference level to be 
the free surface (a/D=l) , vertical suspended sediment 
concentration distributions may be extrapolated as a 
function of the longitudinal variation of suspended 
sediment at the surface. 

F. A prescribed flux is specified at the channel bottom 
to simplify modeling the sorption influence of bedload 
and bottom sediment. As previously mentioned, the 
effects of bedload and bottom sediment on the pollutant 
mass transport process have been investigated in two 
separate studies, hence, an interesting line of research 
would be to integrate the effects of bottom sediment, 
bedload, and suspended sediment into a single convective- 
diffusion model in an attempt to obtain a more detailed 



understaBding of the pollutant mass transport process. 
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LIST OF SYMBOLS 

a = empirical constant equation (2.13) 

ae = artificial diffusion coefficient 

\ 

a^ = constant equation (3.6) 

a 2 = constant equation (3.6) 

a . 3 = constant equation (3.6) 

a 4 = constant equation (3.6) 

as = constant equation (3.6) 

b = empirical constant equation (2.13) 

c = courant number 

C = concentration of pollutant 

= initial pollutant concentration 
Cas = concentration of sorbed pollutant 

Ce = equilibrium concentration of sorbed pollutant 

Co = constant of integration 

Cs = suspended sediment concentration 

D = depth of flow 

= turbulent vertical diffusion coefficient for 
^ suspended sediment 

ds = mean suspended sediment particle size 

E - longitudinal dispersion coefficient 

~ vertical momentum exchange coefficient 

Emag = value of longitudinal diffusion coefficient 
prescribed at middepth 
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= longitudinal diffusion coefficient 
E = nondiinension longitudinal diffusion coefficient 

6y = vertical diffusion coefficient 

Ey = nondimensional vertical diffusion coefficient 

g = acceleration of gravity 

G = energy dissipation rate per unit mass of fluid 

I = subscript denoting longitudinal cartesian direction 
J = subscript denoting vertical cartesian direction 

k = empirical Freundlich constant 

kg = sorption rate constant 

= nondimensional sorption rate constant 
L = characteristic eddy length 
m = empirical Freundlich constant 

M = weight of suspended sediment 

n = Manning's roughness coefficient 

N = subscript denoting time 

r = sorption function 

r* = nondimensional sorption function 

R = hydraulic radius of the channel 

Ri = Richardson's Number 

s = dummy argument equation (3.5) 

So=Se= slope of the channel 
t = time 

T = nondimensional time 
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u = longitudinal velocity component 
U = nondimensional Longitudinal velocity component 
13* = shear velocity 

V = vertical velocity component 

V = average longitudinal velocity 
W = weight of sorbed pollutant 

X = longitudinal cartesian direction 
X = nondimensional longitudinal cartesian direction 
y = vertical cartesian direction 

Y = nondimensional vertical cartesian direction 
= average particle fall velocity 

Z = exponent equation (2.28) 

a = coefficient equation (2.18) 

a' = coefficient equation (2.22) 

3 = coefficient equation (2.20) 

Y = specific weight of water sediment mixture 
Y^ = specific weight of water 

At = time increment 
Ax = longitudinal space increment 
Ay = vertical space increment 

K ~ von Karman's constant 

5 = dummy variable of integration equation (3.5) 

p = density of the fluid 

T = shear stress 

Tq = bottom shear stress 
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Dispersion Experiments 
n = .01 .03 .06 

Eraag = 10 50 100 , 10 50 100 10 50 100 (ftysec) 


Stratification Experiments 


Ca = 100 500 1000 (ppm) 



Sorption Experiments 


Mercury 

DDT 

Heptachlor 


1 X lO'** 

1 X 10‘= 

1 X 10'® 


5 X 10'“ 

1 X 10'“ 

1 X 10'“ 

(1/sec) 

lx 10'^ 

1 X 10'^ 

1 X 10'® 



Figure 4.1 Model Simulations 
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Figure 4.2 Vertical velocit/ distribution as a 
function of Manning's roug-uiess 
coefficient n. 









Figure 4.4 Vertical suspended sedunent concentration 
as a function of tne reference suspended 
sediment concentration Ca. 
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Figure 4.7 Vertical distribution of Sichardson's 

number as a function of the mean particle 
size ds. 
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Figure 5,2b Longitudinal pollutant concentration profile as a function of 
the magnitude of the longitudinal diffusion coefficient 



Figure S,2c 


Longitudinal pollutant concentration profile as a function of 
the magnitude of the longitudinal diffusion coefficient 
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Figure 5.4b Vertical pollutant concentration profiles 
illustrating tbe effect of varying mean 
particles size on tne infl..ence of 
stratification . 
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Figure 5.B I.ongitudinal pollutant concentration profiles illustrating the uptake 
of DDT as a function of the sorption parameter ks. 
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Figure 5.10 Longitudinal pollutant concentration profiles illustrating the uptake 
of lleptaolilor as a function of the sorption parameter ks . 


